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(57) ABSTRACT 

An embodiment of the present invention uses a single 
detection system to approximate a location of lightning 
strikes. This system is triggered by a broadband RF detector 
and measures a time until the arrival of a leading edge of the 
thunder acoustic pulse. This time difference is used to 
determine a slant range R from the detector to the closest 
approach of the lightning. The azimuth and elevation are 
determined by an array of acoustic sensors. The leading edge 
of the thunder waveform is cross-correlated between the 
various acoustic sensors in the array to determine the 
difference in time of arrival, AT. A set of AT’s is used to 
determine the direction of arrival, AZ and EL. The three 
estimated variables (R, AZ, EL) are used to locate a probable 
point of the lightning strike. 

27 Claims^ 8 Drawing Sheets 
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SINGLE STATION SYSTEM AND METHOD 
OF LOCATING LIGHTNING STRIKES 

This application claims the benefit of provisional appli- 
cation No. 60/217,601, filed Jul. 11, 2000. 5 

ORIGIN OF THE INVENTION 

The invention described herein was made in the perfor- 
mance of work under a NASA contract and is subject to the 
provisions of Section 305 of the National Aeronautics and 
Space Act of 1958, as amended, Public Law 85-568 (72 Slat. 
435; 42 U.S.C. §2457). 

TECHNICAL FIELD OF THE INVENTION 

The present invention relates generally to identifying 
locations of lightning strikes. 

BACKGROUND OF THE INVENTION 

Electronic equipment is susceptible to damage caused by 
nearby lightning strikes. The accurate knowledge of a light- 
ning striking point is important to determine which equip- 
ment or system needs to be tested following a lightning 
strike. Existing lightning location systems can provide cov- 
erage of a wide area. For example, a lightning location 
system can provide coverage of an area having a 30 km 
radius. This system, however, has a 50% confidence region 
of about 500 meters. That is, the system has a 50% confi- 
dence that a lighting strike is within 500 meters of an 
identified location. As such, present lightning location sys- 
terns cannot be used to determine whether a lightning strike ^ 
occurred inside or outside of a parameter of an area of 
concern. One such application of a lightning location system 
is a Space Shuttle launch pad for the National Aeronautics 
and Space Administration (NASA). By accurately determin- 
ing lightning strike locations, electronic equipment located ^ 
within the launch pad area can be tested and/or reset to avoid 
erroneous operation. 

One method of determining the location of lightning 
strikes uses a set of video cameras that are pointed in ^ 
different directions within the ai*ea of concern. If a lightning 
strike occurs within the field of view of three or mere 
cameras, the location of the strike can be determined. 
How^ever, if the cameais are not pointed in the correct 
direction, or either an object or a heavy rain downpour 
obscures their field cf view, it is difficult or impossible to 
accurately determine a striking point of the lightning. 
Further, this method has a relatively large uncertainty and 
does not facilitate an accurate location of the exact point of 
contact to the ground. 

For the reasons stated above, and for other reasons stated 
below which will become apparent to those skilled in the art 
upon reading and understanding the present specification, 
there is a need in the art for the system and method to 
accurately approximate locations of lightning strikes. 55 

SUMMARY OF THE INVENTION 

The above-mentioned problems with lightning strike 
location and other problems are addressed by the present 
invention and will be understood by reading and studying 60 
the following specification. 

In one embodiment, a system to determine an approxi- 
mate location of lighting strikes comprises an electric field 
sensor, a plurality of acoustic sensors collocated with the 
electric field sensor, and a processor to determine the loca- 65 
tion of lighting strikes based on a time differential between 
an electric field pulse and a acoustic wave. 
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In another embodiment, a system to determine an 
approximate location of lighting strikes comprises an elec- 
tric field sensor, a plurality of acoustic sound sensors col- 
located with the electric field sensor, a temperature sensor, 
and a processor. The processor determines the location of 
lighting strikes based on a time differential between an 
electric field pulse and a acoustic wave. The processor uses 
a look-up table to compensate for ambient temperature 
measured by the temperature sensor. 

A method for determining an approximate location of a 
lightning strike comprises collecting lightning strike 
information, including a time of arrival of an electric field 
pulse using an electric field sensor and time of arrival data 
of an associated sonic wave from a lightning strike. The 
sonic wave data is detected using a plurality of acoustic 
sensors collocated with the electric field sensor, and pro- 
cessing the lightning strike information to approximate the 
location of the lightning strike. 

Another method for determining an approximate location 
of a lightning strike comprises detecting an electric field 
signal, detecting a sonic pulse with a plurality of collocated 
sensors, measuring a time of arrival difference between the 
electric field signal and the sonic pulse, and determining a 
range, azimuth and elevation of the lighting strike. 

BRIEF DESCRIPTION OF THE DRAWINGS 

FIG. 1 illustrates one embodiment of a single lightning 
location system; 

FIG. 2 illustrates a shape of a wave front distorted by 
wind; 

FIG. 3 illustrates a baseline and hyperbola; 

FIG. 4 illustrates a baseline and hyperbola with an infinite 
source; 

FIG. 5 illustrates a baseline and a location cone; 

FIG. 6 illustrates a definition of AZ and EL; 

FIG. 7 shows a close-up of three intersecting lines of 
position; 

FIG. 8 is a conceptual diagram of an error ellipse; and 

FIG. 9 is a flow' chart of one method of the present 
invention. 

DETAILED DESCRIPTION OF THE 
INVENTION 

In the follovdng detailed description of the prefeited 
embodiments, reference is made to the accompanying draw- 
ings that form a part thereof, and in which is shown by way 
of illustration specific preferred embodiments in which the 
present invention may be practiced. These embodiments are 
described in sufiBcient detail to enable those skilled in the art 
to practice the present invention, and it is to be understood 
that other embodiments may be utilized and that logical, 
mechanical and electrical changes may be made without 
departing from the spirit and scope of the present invention. 
The following detailed description is, therefore, not to be 
taken in a limiting sense, and the scope of the invention is 
defined only by the claims. 

The present invention provides a system and method that 
can accurately locate a lightning strike within an area of 
interest. In one embodiment, the present invention can 
accurately locate a lightning strike within a few meters. The 
fast varying electric current associate with lightning dis- 
charges generates large electric field variations. The electric 
field waveform propagates at the speed of light in a radial 
direction from the striking point of lightning. The sudden 
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heating of the air caused by the large currents associated 
with the lightning discharge produces a sudden expansion of 
the air near the lightning channel. This results in a sonic 
wave (thunder) that initially, for the first few meters, propa- 
gates at a supersonic speed and later propagates at a sonic 5 
speed. 

For an observer located remotely from a lightning strike 
location, the electric field waveform arrives earlier than the 
sonic sounds. This is because the electric field waveform 
travels and a speed of approximately 300,000, (X)0 m/s, while jq 
the sonic wave travels and approximately 350 m/s. An 
observer can determine the distance to the striking point by 
measuring the time difference between the arrival of the 
electric field waveform and the arrival of the sonic wave. 
Tliis measurement defines a circle, with the observer at the 
center, in which the lightning strike might have occuned. A 
second observer at a different location using the same type 
of measurement also has a circle defined around them in 
which the lightning might have occurred. These two circles 
intersect at two points. With the addition of a third observer, 
a single striking point can be determined. 20 

In an embodiment described in U.S. patent application 
Ser. No. 09/784,405 filed Feb. 13, 2001 now U.S. Pat. No. 
6,420,862, a lightning locator uses several (at least three) 
remote receiver stations each equipped with a broad band 
antenna and a microphone. The antenna detects lightning 25 
through its strong electromagnetic pulse and starts a timer. 
Tlie timer is turned off by a detection of the arrival of the 
sound wave due to thunder. Range is estimated by multi- 
plying the time delay by the speed of sound. The ranges 
d(5tected by each receiver are transmitted to a central pro- 3Q 
ce^ssor to estimate the location of the lightning strike. 

An embodiment of the present invention uses a single 
detection system. This system is triggered by a broadband 
RF detector and measures a time until the arrival of a leading 
edge of the thunder acoustic pulse. This time difference is 35 
used to determine a slant range R from the detector to the 
closest approach of the lightning channel. The azimuth and 
elevation (AZ and EL) are determined by an array of 
acoustic sensors with one to several meter baselines. The 
leading edge of the thunder waveform is cross-correlated 40 
b<;tween the various microphones in the array to determine 
the difference in time of arrival (DTOA or AT). A set of AT’s 
is used to determine the direction of arrival (AZ and EL). 
The three estimated variables (R, AZ, EL) are used to locate 
a probable point of the lightning strike. 45 

The present invention provides a lightning location sys- 
tem that uses a single station to determine an approximate 
location of lightning strikes. One embodiment of the present 
Icjcation station 100 includes an electric field sensor 102 and 
five sonic all receivers 104fl-e, see FIG. 1. Four sonic 50 
sensors 104a-i/ are located in a periphery of a circle at 90 
degrees from each other. For example, the sensors can be 
placed on a 5-meter radius circle. The electric field sensor 
102 is located at the center of the circle. The fifth sonic 
sensor 104e is located above the plane formed by the four 55 
s(msors, such as 2 meters above the plane. The station, 
therefore, resembles a pyramid or a conical shape. A pro- 
cessor 106 is coupled to the sonic and electric field sensors 
to control the station and process data. As explained herein, 
mnd speed sensor 108, temperature sensor 110 and a 60 
humidity sensor 112 can be provided to correct the data 
processing and reduce errors. A camera 114 can also be 
added to provide a visual reference during operation. The 
output from the system can provide both data output and 
visual output. 65 

In operation, the electric field sensor receives an electric 
field within a few microseconds after a lighting strike. A 
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sonic wave arrives at the sonic sensors following the electric 
field. The sonic field reaches each sensor within millisec- 
onds of each other, since the distance to the striking point is 
different for each sonic sensor. Based on the differences in 
the time of arrival of the sonic wave at each sensor, a set of 
equations is solved to determine the angle of arrival of the 
sonic waveform The differences in the time of arrival can be 
precisely measured by performing digit cross-correlations 
among the signals received by the five sensors. The distance 
to the lighting strike is determined by the time elapsed from 
the detection of the electric field to the detection of the 
thunder signal (sonic). The speed of sound can be estimated 
by monitoring the ambient temperature and automatically 
compensating for any deviations during the course of the 
day. A look-up table is used to determine the correct sound 
speed to use for the calculation of the distance to the 
lightning strike. 

Digital signal processing techniques are used to discrimi- 
nate between thunder waveforms from a close lightning 
strike (within about a 1-mile radius) and thunder waveforms 
from a preceding distant strike. This process is important 
since, during the active phase of a thunderstorm, lightning 
activity could occur at a rate exceeding 10 strikes per 
minute. The discrimination between close and distant thun- 
der is based on the frequency content of the thunder wave- 
form. Propagation through the air causes the high-frequency 
components to attenuate much faster than the lower frequen- 
cies. Thunder generated within short distances has typical 
high frequency components (a “clap” of initial sound) 
followed by the longer duration rumble. Thunder from more 
distant waveforms lose their high-frequency components 
and no “clap” sound is heard. 

A designer of a single station sonic lightning location 
system can select the geometry of the sonic sensor micro- 
phones in a variety of manners depending on the portion of 
the sky that is of greatest interest. Ttie analysis provided 
herein allows a designer to calculate the accuracy of direc- 
tion of arrival estimates based on user provided locations of 
the sensor microphones and the timing accuracy possible. 
The timing accuracy is based on the electronic and digital 
design of the detection system and also on the nature of the 
thunder waveforms actually experienced under field condi- 
tions. The math presented provides the estimated azimuth 
and elevation of the source and estimates of the confidence 
levels of these estimates (enor ellipse). 

The practically instantaneous arrival of the electromag- 
netic pulse is detected with the broad band antenna 102 
which starts a timer in the processor 106. When the thunder 
signal arrives at a selected microphone 104, sonic sensor, the 
timer is stopped. A variety of techniques for digital signal 
processing can be used to start and stop the timer and will 
not be addressed here. A signal analysis technique is used to 
select only nearby events by filtering, see U.S. application 
Ser. No. 09/784,405 filed Feb. 13, 2001. Any known tech- 
nique can be used to estimate a time delay of arrival of a 
sonic or thunder pulse accurate to within a few microsec- 
onds. This time delay is multiplied by c, the speed of sound, 
to estimate the slant range to the lightning attach point 
(“ground” contact). With the speed of sound of roughly 300 
meters per second, this provides estimation accuracy on the 
order of a centimeter or less. However, the variation of the 
speed of sound with temperature and humidity and the 
effects of wind degrade this ideal performance. 

Wind velocities are significantly lower than sonic veloc- 
ity. For example a wind speed of 50 miles per hour equates 
to 22.6 meters per second, only 7.5% of the speed of sound 
(approximately 300 m/s). The thunder leaves the lightning 
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channel as a cylindrical wave front with a velocity that is the 
vector sum of the wind speed and the sonic speed. This 
affects the overall travel time to the sensor, used to estimate 
the slant range, and potentially, the curvature of the wave 
front at the sensor affects the estimated azimuth and eleva- 
tion (for the single station system). Assuming that the 
lightning channel is vertical and the wind velocity vector is 
strictly horizontal, any effects to the estimated elevation 
angle can be ignored. This reduces the problem to two 
dimensions. 

As shown in FIG. 2, a lightning location point is given by 
P and the sensor location is at the origin O. The wind 
velocity (to the south in the figure) distorts the wave front by 
adding the components of the wind velocity to each wave 
front velocity vector. In the presence of wind, the velocity 
vector of the wave front measured along an azimuth 0 is 
given by: 

V (0)='?-| V wlcos(0-3) ( 1 ) 

The second term of equation (1) measures the contribu- 
tion of the wind along the line of sight to the event P. The 
minus sign comes from consideration of the angles involved. 
Note that when 3 is greater than 0 (as in FIG. 2) the sign of 
the cosine term changes to plus increasing the effective 
wavefront velocity. Thus, in the presence of wind, the wave 
front arrives sooner or later than it would in still air. 
Designating “AT’ as the time delay between the lightning 
flash and the thunder arrival and recognizing that wind has 
no effect on the electromagnetic propagation time, the range 
is: 


difference is the difference time of arrival (DTOA), in time 
units, times the speed of sound c. The location of S is 
unknown and, hence, the orientation of the plane, S lies on 
a “hyperbola of revolution” or a hyperbolic cone given by 
5 rotating the hyperbola around the axis defined by the base- 
line. 

The curve of a hyperbola, away from the origin, asymp- 
totically approaches a straight line that deviates from the 
baseline with interior angle \)/, see FIG. 4. In the present 
10 invention, only distant points relative to the length of the 
baselines are used. Therefore, the equations are simplified 
by using only straight lines instead of hyperbolic curves. 
When the straight line is rotated about the axis, it generates 
a circular cone (FIG. 5). 

*5 The present system measures a time difference in micro- 
seconds between the time of signal arrival at one micro- 
phone minus the time of arrival to a second microphone. 
This DTOA times c, the speed of sound, is a distance 
difference. The measured TOA from microphone 1 is Tl, 
20 and the measured TOA from microphone 0 is TO then: 

Aro,=7D-n (5) 

If the ATqi term is positive, then TO is a larger number and 
25 is later in time; the event reached “1” first. Therefore, the 
cone angles toward microphone 1. Likewise, if ATq^ is 
negative, the event reached “0” first and the cone opening 
lies toward microphone 0 (\|c>7c/2). The cone of constant 
ATqi forms a constant angle between the baseline and the 
30 cone surface y (FIG. 2). In terms of the ATq,, for sources far 
away from the baseline: 


R^^T V^-^T{c-Vw cos(O-a)) (2) 

This equation must be inverted so that the measurement 
AT is expressed in terms of unknown variables: 


cAroi — /qi cos y 

35 So that the cone angle y is given by: 


( 6 ) 


Ar = 


R 

(c- Vwcos(6-S)) 


( 3 ) 


^ = arc cos 




( 7 ) 


The value of c changes with temperature and humidity 
(molecular weight of the atmosphere). The value of sound 
speed in a gas is given by: 

c=iyRor/My’'^ ( 4 ) 

45 

where y is the ratio of specific heats, Rq is the universal gas 
constant (8314 J kg“* K"'), T is the absolute temperature and 
M is the molecular weight. The speed of sound varies as the 
square root of the temperature. The speed of sound increases 
by 0.6 m/s for each degree increase in Celsius temperature. 50 
Tiis significant effect should be corrected using temperature 
sensor 110. The presence of water vapor causes both y and 
M to decrease. Tie water vapor effect can often be ignored 
since its correction is typically less than 1.5%. However, if 
humidity data is available to the system it would be advis- 55 
able to correct for its effects. To preserve accuracy, an 
external measurement of local wind speed and temperature 
may be desired. Measuring humidity provides a smaller 
benefit. If such measurements are provided the range accu- 
racy of the present single station system is quite high. 60 
Any two sonic receivers (R) form a “baseline,” that is a 
line connecting the two receiver coordinates. The location of 
a source S that produces a signal lies on a hyperbola in the 
plane formed by two microphones and the lightning source, 
see FIG. 3. For each point on a hyperbola, the difference 65 
between the distances to the two foci is a constant. Eacii 
microphone forms a focus of the hyperbola. This distance 


This cone traces a circular arc on the sky. Using multiple 
baselines can be used to produce multiple arcs, such ihat 
their intersection defines the azimuth and elevation of the 
source. In order to accomplish this, the baseline orientation 
is defined in terms of azimuth and elevation, which, coupled 
with the above DTOA equation, allow the definition of the 
arc. Tie AZ and EL angles can be defined in terms of local 
Cartesian coordinate system. The orientation of each 
receiver pair baseline and the direction to each source in 
terms of AZ and EL is expressed as follows. The locally flat 
Cartesian coordinate system is called (east, north, up) for (x, 
y, z) respectively. The azimuth AZ and elevation EL angles 
are defined with reference to FIG. 5. Assuming that the 
reference receiver (point 0) lies at the point (xO ,y0, z0)=(0, 
0, 0). The baseline vector is defined as the vector from (0, 0, 
0) to the receiver listed first in the DTOA equation. So that 
the baseline “10” vector runs from the point (0, 0, 0) to (xl, 
y 1, zl). This vector defines a baseline orientation in terms of 
azimuth and elevation. The azimuth is measured from the 
north (y axis) in a clockwise direction. The elevation angle 
is measured up from the (x, y) plane. Points lying in the (x, 
y) plane define position vectors with a zero elevation angle. 

The baseline vector is therefore given by: 


^ io=(Jf i-0)/+(yi-0) 

}Hzi-0)k 




f 


( 8 ) 
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The length of this vector is given by: 


11/ 10 II = + A + =/ io . 


(9) 


where Ijq, I 20 I 30 indicate the length of three baselines 

vectors. It must be kept in mind that the results below will 
depend on the sign of the DTOA, which is implicitly affected 
b> the direction assumed for the baseline vector. Given a 
point (x^, y^, z^), which could be an antenna or a source 
point, it is transformed to spherical coordinates as follows, 
where again 1 denotes the length of the vector from the origin 
to the point, see FIG. 5: 


10 


Zo=l sin EL 


15 


Xo=l sin AZ^l cos EL sin AZ 
yo=l cos AZ=l cos EL co&AZ 


( 10 ) 


Since we are only interested in AZ, EL, we divide through 
b)' the length and solve for the angles; 


20 


£L = arcsin^^j 


(11a) 


25 


AZ = arc t 




(lib) 


The angles of the baseline vector can be defined using 
receiver 0 as the reference origin ( 0 , 0 , 0 ): 


30 


ELio = arc sin 


■ f 
in 


35 


AZio = arctan|^ 


(12) 


40 

Equation (12) defines the AZ and EL of any point in terms 
of our Cartesian system. 

'Fhe following section derives the equation for the arc on 
the sky based on a measured delta T, baseline length, and 
baseline orientation, and is a general treatment which can be 45 
solved for any baseline orientation. Given a source point 
S(x^, y^, z^) or, in spherical polar coordinates, ( 1 „ AZ^, El^) 
and a baseline unit vector; 


(^(TO-r/)) 


(13) 50 


Expressing the vector pointing from the origin to the 
source point: 35 


s=xj+yj+zjc 

and its corresponding unit vector; 


r, r, r. 


(14) 


(15) 


60 


Note that the dot product between these two vectors is 
equal to the cosine of the angle between them. 65 
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The dot product is given by the sum of the products of the 
three components; 


5../ = , Zl + . £i 

rs ri rs n rs n 


(17) 


Using the following spherical coordinate conversion for- 
mulas (equation 10 ) to express the dot product in terms of 
angles; 

s l=cos EL^ sin AZ, cos EL^ sin AZ, +cos cos AZ, cos £L, cos 

AZ,+sin £L, sin ££, (18) 

Rearranging; 

S‘l=cos ELj. cos ££, sin AZ, sin AZ, +cos £L, cos ££, cos AZ, cos 
AZ,+sin £L, sin ££, (19) 

Using 

cos^ = (^(TO-ry)), (equation?) 


the equations for an arbitrary baseline becomes: 

-^(ro - Tl) = cos ££5 Cos ££isin AZs sin AZi + 

/oi 

cos ELscos £Licos AZ^cos AZi + 
sin££5sin££i 


Each baseline in the system produces an equation 20. 
Solving the set of equations for an array of baselines is 
addressed below. In order to solve for the azimuth and 
elevation of a distant source, at least three baselines are 
needed. The lines of position provided by each baseline 
describe latitude type lines (i.e. not great circles) across the 
sphere of the sky. Two baseline-derived lines of position 
(LOP’s) intersect at two points. A third baseline is required 
to resolve the ambiguity. 

Four sensors are needed to obtain three baselines. Care 
must be taken in selecting the baselines produced by the four 
sensors. In particular, three independent sets of baselines 
must be chosen such that one baseline cannot be expressed 
as the difference between two others. This causes the matrix 
equations to have a lower rank than three and be under 
determined. This effect prevents production of three base- 
lines from three sensors. Beyond this constraint, however, 
the system can actually select the arrangement of baselines 
to provide the best solution for each particular event, if the 
system geometry allows it. Therefore, all of the solutions are 
over-determined in the sense that there are more equations 
than unknowns. There may also be constraints such as the 
ground level that can restrict solutions to portions of the sky 
(so that elevation angles can only be between zero and 
positive pi over two). 

Since each baseline measurement contains errors, it is 
desired to find the solution AZ and EL point that minimizes 
the sum of the square errors. The solution point is called 
(AZ£ 57 > ELg 57 ^) and using this point the expected AT’s for 
each baseline can be computed; 


AT) = -(cos£L,cos£L,sinAZ,sinAZ; + 
c 

cos ££jCOS ELiCos AZ^cos AZt + 
sin £L,sin £D 


(17) 


The angles that involve the i** baselines can be defined as 
constants as follows; 

a/=cos ELi sin AZ„ p,=cos ELf cos AZ„ ypsin £L 


i'l=cos Y 


(16) 


(18) 
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With cos EL equal to CEL and sin AZ is SAZ, these 
changes simplify equation 17 to: 

^Ti = -{OiCELSAZ-^ ^iCELCAZ + nSEL) 
c 

A nonlinear regression on this set of equations produces 
a least squares estimate of the (AZ, EL) of the source. An 
initial very rough estimate of the source location must first 
be produced. This can be done in several ways but a method 
is recommended whereby the closest grouping of three 
intersections is identified. This can be accomplished by 
creating a table of computed delta T*s for each 1() degree by 
10-degree patch of sky. Since AZ runs from 0 to 360 degree 
and EL runs from 0 to 90 degrees, this results in a listing of 
324 delta T sets (3 values each). The actual delta T’s can 
then be subtracted from each computed delta T then squared. 
By finding the smallest value in this list, an estimate is 
provided of the correct value. Using this as a starting value, 
the delta T functions can be extended in a Taylor series 
around the estimating starting point. By using only the first 
order terms of the series, a set of linear equations is obtained 
as follows: 


^Eueasuredi ~ ATEstimatedi + 


dinTi 

dAZEst 


{AAZ)-h 


dATj 

dELEsf 


(AEL) 


( 20 ) 


25 


The partial derivative terms are computed as follows: 

30 


OAZes, 


U 


{QiCELE^ CAZest -fiiCELEsr SAZe^) 


(21a) 


dTi 

dELEst 


—{-QiSELEst SAZesi - PiSELesi CAZesi + 


yiCELEst) 


35 

(21b) 


Since the values of all of the parameters are known, these 
terms reduce to numbers for any particular choice of esti- 
mated AZ and EL position. With 3 (or more) baselines, a 
matrix, equation is coastructed which is linear in the differ- 
ences in the delta T"s and the corrections to AZ and EL. 
With: 


^AT r AT Estimated i 

the matrix equation can be written as: 


'SATi ^ 
8 AT 2 
.SATi, 


dTi 

dT, 

dAZEs, 

dELEst 

dTj 

dT2 

dAZEst 

dELEst 

dTi 

dTi 

dAZEst 

dELEst 


( 22 ) 
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(23) 
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This can be written in matrix vector notation as: 

SAT=HAe (24) ^ 

This is solved using the generalized inverse (multiply 
both sides by the transpose of matrix H then take the inverse 
of the now square H^H). The solution is written as follows: 

AB=(H^Hr'H^^T (25) 

o5 

The values of delta AZ and delta EL are added to the 
original estimated values of AZ and EL to produce a new 


estimate. The values of the squared 6AT terms are summed. 
The process is repeated with the new estimate until the sum 
of the squared residuals no longer decreases. This represents 
an un- weighted least squares estimate. 

An optimal solution would weight each baseline based on 
the inherent accuracy provided. Returning to the original 
equation for the delta T as a function of cosine of the angle 
between the baseline and the direction to the source: 


and taking the derivative of both sides (leaving 1 constant of 
course): 


dAT = --sini/zdiff 

c 


(27) 


Rearranging to find the error in angle (d^) with error in 
delta: 


dt/f ^ (28) 

dAT Isioif/ 


This error function is a minimum when \|/ is at 90 degrees, 
when the source is directly broadside to the baseline. The 
function has a singularity at X|/ equal to 0 or 180 degrees. A 
weighing function is selected that gives more weight to 
baseline measurements taken when delta T is very small and 
giving little weight to measurements where delta T is near its 
maximum value (1/c). This relationship is expressed in terms 
of variances by squaring: 


al = 




(29) 


By using a trigonometric identity, the angle y is elimi- 
nated and delta T is used instead: 



When combining measurements from various baselines, it 
is desired to weight the results to favor those baselines with 
lesser errors. To do this, a weight matrix W is introduced. 
This matrix has the properties that it is diagonal (which 
assumes measurements are uncorrelated) and that its trace is 
equal to one (does not change the value of the final result). 
Assuming that the basic errors in timing determination are 
the same for each baseline, a fixed value to the timing 
variance is assigned. The measurements can be weighed as 
follows: 


( 1 




(31) 


This weight matrix is employed to obtain a weighted least 
square estimate by multiplying both sides of the model 
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equation then rearranging into the generalized inverse. This 
proceeds as follows: 

WbAT^WHAe (32) 

Solving for the angle adjustments: ^ 

AO=(H^WH)-^H^WbAT (33) 

Using this equation, the measurement variances can be 
propagated into the covariance matrix of the AZ, EL angles lo 
solution. The covariance matrix of measurements is com- 
posed of the variances of the measurements on the diagonal 
and the covariances or correlations off diagonal as follows: 


systems and perform studies regarding accuracies required 
and so on. This also allows the computation of the expected 
accuracy for specific detected events and to describe the 
error ellipse onto an image of the scene to assist in identi- 
fying the location of the lightning attachment. 

A question arises regarding the change in apparent direc- 
tion to the source caused by the wind. This effect could be 
serious in the case of the present single station system. An 
embodiment of the sensor array is made of individual 
detectors spaced a distance D apart, where D is a charac- 
teristic dimension rather than an exact dimension (for esti- 
mation purposes). At a range R, a distance D subtends an 
angle d0. Since d0 is very small: 

D^dQ (42) 


COV(ATi)^E{EAr.elr,] = I>£,T = 


CTj (T 12 fTi3 ' 
^21 

CT3I (T32 CT^ ^ 


(34) 


Since it is assumed that the errors are the same for each 
baseline and they are not correlated, this reduces to the 
identity matrix times a single measurement covariance: 


15 The following analysis shows that the wind does not 
significantly affect the measured angle of arrival. If the 
velocity as a function of 0 is suflScient to cause an additional 
difference in time of arrival between the two sensors then the 
derived angle to the source will be incorrect. By taking the 
20 derivative of the velocity expression as a function of angle 
0, the time difference measured can be derived. The deriva- 
tive of V with angle is: 




(35) 


The covariance matrix of the angle solution is found by 
using the weighted generalized inverse equation above: 


25 




(36) 


The covariance matrix for the angles is given by the 30 
following: 


dV (43) 

— =Vwsvn{G~d) 


With an expression for d0 in terms of R and D: 


DVw (44) 

AV= — ^sin(0-(» 




^AZ 


O' AZEL 


O'ELAZ ^ 


Applying this expression to the expression involving AT 
and velocity, the additional AT or the 5AT as a function of 
35 array dimension, range and wind velocity can be estimated 
as: 


This matrix is positive symmetric and therefore defines an 
ellipse on the sky centered around the estimated source 
point. The ellipse has semi major and semi minor axes and 
inclination to the azimuth axis are defined by: 



+ ^AZEL 


^(^^AZ ” <^Elf + <^AZEL 


tan2y = 


'^AZEL 


(39) 


(40) 

(41) 
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6AT = -^AV^ 


g-Vw 
c^sin {9-d) 


(45) 


With an array dimension of one meter, a wind speed of 50 
mph (22.6 m/s), a range of 100 meters, and the sine function 
equal to one (wind along direction to source); a timing skew 
of only 0.00025 seconds is found. This can be translated into 
angular error by remembenng that for a baseline D between 
two receivers that measure a time of arrival difference of AT, 
the angular error caused by an error in AT is: 


\DlsiTut 


(46) 


Note that the angular portion of the sky is being treated as 
being locally Cartesian which is valid since this is a small 
patch of sky. If the input values for timing error represent 
one standard deviation of a random sampling of such timing 
eirrors then the sigma’s above also represent a similar 
confidence region (i.e. 68.3% of expected events occur with 
±G of the estimated point. However, the ellipse, because it 
njpresents cross correlations will enclose only 39.4% of 
events for the same 1 sigma input. This low confidence level 
may cause some users to utilize an ellipse enclosing a 95% 
confidence limit, which requires use of 2.447 times one 
sigma for the measurement input. 

The tools presented herein allow one to compute the error 
ellipse properties for all portions of the sky given the 
configuration of the baselines (lengths, angles, timing errors) 
and thereby compare different configurations of detection 


The angle y is used to differentiate this angle from the 
azimuth to the source (which is derived from a number of 
baseline angles). Inserting in the value for \|/ as 45 degrees 
(sine equals 0.7071) and using the previous numbers for 
timing error and so on, an angle error of 0.106 radians or 
0.0019 degrees or 6.7 seconds of arc is determined. 

In the above example, the range error was roughly equal 
to the ratio of the wind speed to the sound speed which is 
significant. For the case of angular determination, the error 
was so small that it can be neglected. This means that the 
single station system may require and external source of 
wind speed and direction to correct the measured range 
values (or be accepted as a known systematic error). 

In an alternate method, using an array of three baselines 
of equal length and mutually orthogonal (such that we have 


55 


60 


65 


US 6,552,521 B1 


13 

a baseline along our x, y, and z axes), an additional equation 
of constraint can be provided. The following sections indi- 
cate potential uses for that equation in assisting in the 
solution of direction finding problems. Equation 2 translates 
the AT measurements into the direction cosines between 
each axis and the vector pointing to the source. Remember- 
ing that the direction cosines for an orthogonal array follow 
the relationship: 


a+cos^ Y 


(47) 


where a, p, and y are the angles between the vector pointing 
to the source s and the x, y, and z baselines respectively. 
Equation 47 is used to perform a quick validation of the raw 
AT measurements, which is a key advantage of an orthogo- 
nal array. For orthogonal baselines: 


Baseline 

Cos ELq, 

Sin EL), 

Cos AZq, 

Sin AZoi 

01 (X) 

1 

0 

0 

1 

02 (y) 

1 

0 

1 

0 

03 (z) 

0 

1 

0 

0 


Using these values in equation 17 results in three mea- 
surement equations relating the AT’s to AZ and EL: 


(7'o-7’,)=cos EL sin AZ 
(Tq-T2)=cos el cos AZ 
{TQ—T-^)=sm EL 


^(^T\/£,45--Ar7^^ouy)^Mininiunri 


- CELSAT)^ + (A7b2 - CELCAZf + (AT’o., - SEL)^] = 0 

oAZ 


(50) 
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This is because ATq 3 gives a line of position across the sky 
that is at a fixed EL and has no dependence on AZ. 

Each measurement provides an arc that makes a large 
circle on the sky. If the measurements contain some errors, 
5 the three circles will not intersect at a single point but will 
form a triangle as shown in FIG. 7. Each comer of the 
triangle represents the AZ and EL solution found using only 
two of the three baselines. The question is where in AZ and 
EL is the “best” or most probable solution given the mea- 
10 sured ATs? 

Equations (41) and (42) provide the least squared error 
(i.e., the sum of the distances between the chosen AZ, EL 
point and the lines of position is a minimum). Unfortunately, 
this procedure is influenced strongly by a measurement that 
15 is in error so is not optimum. 

Measurement errors consist of random and bias errors. 
Random errors are due to noise and other random processes 
in the measurement chain. Random errors are typically 
minimized by averaging a number of measurements. This is 
20 not possible in a lightning measurement system since the 
event only occurs once. Bias errors can, at least in theory, be 
eliminated through a process of calibration and alfect, in 
some consistent way, all measurements. Bias errors that are 
unknown are generally treated as random since they are 
25 unknown. 

To find how each measurement may be affected by 
random errors, equation ( 6 ) is examined and the derivatives 
of both sides are taken: 


(48) 


30 


Note that it is assumed that all of the Ts are equal and 
l/c=l . The subscript “s” on the angles is also dropped. At the 
end of the derivations, the 1/c coefficient is reinserted. To 
further simplify, “cos EL” is written as “CEL”, etc. To solve 
(48), a minimum least-squared-error approach is used to find 
AZ and EL so that: 


(49) 


Each of the terms in parenthesis are call “residuals.” To 
find the solution for AZ, a partial derivative of (49) is taken 
with respect to AZ and set to zero: 


cdET=~l sin 
Rearranging: 


di/f _/c^ -1 
dET 


(54) 


(55) 
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This expression can be interpreted as the ratio of standard 
deviations. Taking the absolute value (since the timing error 
can be positive or negative): 


40 


iT^ 


= (7)i 


(T{^T 


(56) 


If the standard deviation of the timing difference error is 
1 |bis, the baseline length is 2 meters, and using a value of c 
equal to 330 m/s, we see that: 


Evaluating this expression: 


(51) 


and evaluating as above: 

AFo3 


EL = tap 


(AToiSAZ+^ToiCAZf 


50 


Note that the arc tangent has ambiguities due to multiple 
values of AZ with the same arc tangent. Solving for EL: 


^ U AToi - CELSAZi^ + (ATtp - CELCAZ)^ + (^To^ - SELf] = 0 
oEL 


(53) 


The equations contain a couple of interesting results. 
First, equation (51) for AZ contains no dependence on AT 03 . 


(0.000165) (\S0\ 
— x| 
sin^ 


■(”) 


(57) 


Where the last term converts from radians to degrees. The 
best accuracy is achieved when the source is at right angles 
to the baseline and is equal to 0.009 degrees (in the absence 
55 of other errors). Note that this is roughly Vi an arc minute. 

FIG. 8 shows the angular accuracy as a ftincdon of \|/ 
(“delta psi denotes a\|f the standard deviation of \|/). This 
indicates that, in areas where a given baseline is parallel to 
the incoming rays, it should probably be excluded from the 
60 solution. 

One embodiment of a data processing algorithm is pro- 
vided below, with reference to FIG. 9, which results in an 
estimate of the range, azimuth, and elevation of the source. 
The system is placed in a ready state 200 and waits for an 
65 electric field trigger event 202, A timer is started 204 upon 
trigger and the system waits for a sonic wave to be detected 
206, If sonic data is not received by a predetermined 
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threshold time, such as 150 fisec, the timer is stopped 212 
and the system reset 230. If sonic data is received 206, the 
data is filtered for distant strikes 208 and the timer is stopped 
210. A range for the strike is calculated 214 and the sonic 
data is cross-correlated 216. The range and sonic calcula- 
tions can be adjusted for environmental conditions 218, such 
as wind, temperature and humidity. The data is validated 220 
b>' taking the sum of the squares equal to “1” (equation 17). 
If the data is not valid, the system is reset. All AT’ s are 
validated 221 that they fall below “0.9.” If valid, AZ and EL 
are calculated 222 using equations 21 and 22. If not valid, 
high AT measurements are excluded and the AZ and EL are 
calculated 224 using remaining equations in (18). AZ and EL 
errors are calculated 226 by setting values of AT higher and 
lower by 1 sigma, and AZ and EL are recomputed. AZ±AAZ, 
EL+AEL, R±AR are output 228 and the system is reset 230 
for a next event. 

The above algorithm is one embodiment of a methodol- 
ogy to approximate the location of a lightning strike, other 
methods can be implemented and the invention is not limited 
to the above methodology. Although specific embodiments 
have been illustrated and described herein, it will be appre- 
ciated by those of ordinary skill in the art that any 
arrangement, which is calculated to achieve the same 
purpose, may be substituted for the specific embodiment 
shown. This application is intended to cover any adaptations 
or variations of the present invention. Therefore, it is mani- 
festly intended that this invention be limited only by the 
claims and the equivalents thereof. 

What is claimed is: 

1. A system to determine an approximate location of 
lightning strikes comprising: 

a single electric field sensor; 

a plurality of acoustic sensors approximately collocated 
with the electric field sensor; and 

a processor to determine a distance to lightning strikes 
based on a time differential between an electric field 
pulse and an acoustic wave at one of the plurality of 
acoustic sensors, and an angle of arrival based on a 
difference of time of arrival of the sonic wave at the 
plurality of acoustic sensors, wherein the processor 
outputs a range, azimuth and elevation to approximate 
the location of the lightning strike relative to the 
location of the system. 

2. The system of claim 1 wherein the processor compen- 
sates for a speed of sound based on a measured temperature 
and/or humidity. 

3. The system of claim 1 wherein the processor compen- 
sates for the time of arrival of the sonic wave using a 
measured wind speed. 

4. The system of claim 1 wherein the processor provides 
an estimate of the accuracy of the approximate location. 

5. The system of claim 1 wherein the plurality of acoustic 
sensors comprise five acoustic sensors. 

6. The system of claim 5 wherein four of the five acoustic 
sensors are located on a perimeter of a circle on a horizontal 
plane, the electric field sensor is located inside the circle, and 
a remaining one of the five acoustic sensors is located 
vtjrtically above the circle. 

7. The system of claim 5 wherein four of the five acoustic 
sensors are located on a common horizontal plane. 

8. The system of claim 7 wherein the electric field sensor 
is located on the common plane and between the four 
acoustic sensors. 

9. The system of claim 8 wherein a remaining one of the 
five acoustic sensors is located vertically above the common 
plane. 
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10. The system of claim 6 wherein the circle has a radius 
of about 5 meters and the one acoustic sound sensor is 
located about two meters above the electric field sensor. 

11. A system to determine an approximate location of 

5 lighting strikes comprising: 

an electric field sensor; 

a plurality of acoustic sound sensors collocated with the 
electric field sensor; 

a temperature sensor; and 

10 a processor to determine a distance to lighting strikes 
based on a time differential between an electric field 
pulse and an acoustic wave at one of the plurality of 
acoustic sensors, and an angle of arrival based on a 
difference of a time of arrival of the sonic wave at the 
plurality of acoustic sensors, the processor uses a 
look-up table to compensate for a speed of sound based 
on ambient temperature measured by the temperature 
sensor. 

12. The system of claim 11 wherein the plurality of 
acoustic sound sensors comprise five acoustic sensors, 
wherein four of the five acoustic sound sensors are located 
on a perimeter of a circle on a horizontal plane, the electric 
field sensor is located above the circle, and a remaining one 
of the five acoustic sound sensors is located vertically above 
the circle. 

25 13. The system of claim 11 wherein the processor further 

compensates for the speed of sound based on a measured 
humidity. 

14. The system of claim 11 wherein the processor com- 
pensates for the time of arrival of the sonic wave using a 

30 measured wind speed. 

15. The system of claim 11 wherein the processor outputs 
a range, azimuth, and elevation relative to the system, 

16. The system of claim 15 wherein the processor further 
outputs an estimate of the accuracy of the range, azimuth, 
and elevation. 

17. A method of determining an approximate location of 
a lightrdng strike comprising the following steps: 

measuring the ambient temperature; 

collecting lightning strike information, including a time of 

4 Q arrival of an electric field pulse using an electric field 
sensor and a time of arrival data of an associated sonic 
wave from a lightning strike, wherein the sonic wave 
data is detected using a plurality of acoustic sensors 
collocated with the electric field sensor; and 

45 processing the lightning strike information including 
compensating for a speed of sound based on the ambi- 
ent temperature approximate the location of the light- 
ning strike. 

18. The method of claim 17 wherein compensating com- 

50 prises accessing a look-up table. 

19. The method of claim 17 wherein processing the 
lightning strike information comprises calculating an esti- 
mate of the accuracy of the approximate location of the 
lightning strike. 

55 20. The method of claim 17 further comprises: 

measuring a wind speed; 

compensating for the time of arrival of the acoustic wave 
based on measured wind speed while processing the 
lightning strike information. 

60 21. The method of claim 17 wherein a processor outputs 

a range, azimuth, and elevation relative to the sensors. 

22. The method of claim 21 wherein the processor further 
outputs an estimate of the accuracy of the range, azimuth, 
and elevation. 

65 23. The method of claim 22 further comprises indicating 

an estimate of the accuracy of the range, azimuth, and 
elevation. 
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24. A method for determining an approximate location of 
a lightning strike comprising the following steps: 

detecting an electric field signal created by the lightning 
strike with a single electric field sensor; 

detecting a sonic pulse created by the lightning strike with 
a plurality of approximately collocated acoustic sen- 
sors; 

measuring a time of arrival difference between the electric 
field signal and the sonic pulse at one of the plurality of 
approximately collocated sensors to determine an 
approximate distance from the lightning strike to the 
electric field sensor; 

measuring a time of arrival difference between the sonic 
pulse at each of the plurality of approximately collo- 15 
cated acoustic sensors to determine an angle of arrival 
of the sonic pulse; and 

providing a range, azimuth and elevation of the location 
of the lightning strike relative to the electric field 
sensor. 

25. The method of claim 24 further comprises compen- 
sating for a speed of sound based on a measured temperature 
and/or humidity. 
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26. The method of claim 24 further comprises compen- 
sating for the time of arrival of the sonic pulse based on a 
measured temperature and/or humidity. 

27. A method of designing a single station lightning 
5 location system comprising; 

positioning a single electric field sensor; 
defining at least a first geometry for locating a plurality of 
acoustic sensors about the single electric field sensor; 
calculating an accuracy of direction of arrival estimates 
based on the at least first selected geometry; 
defining a second geometry for locating the plurality of 
acoustic sensors about the single electric field sensor; 
calculating an accuracy of direction of arrival estimates 
based on the second selected geometry; and 
comparing the accuracy of direction and arrival estimates 
based on the first and second geometries with the 
expected location of the lightning strike to determine 
the geometry providing the greatest accuracy in deter- 
mining the range, azimuth and elevation of the light- 
20 ning strike as compared to the location of the electric 
field sensor. 


